﻿using System;
using System.Text;
using System.Drawing;
using System.Buffers;
using System.Collections;
using System.Collections.Generic;
using System.Runtime.InteropServices;

public static partial class NativeAOT
{
    [UnmanagedCallersOnly(EntryPoint = "mtcl")]
    public static unsafe int mtcl(IntPtr zr_ptr, IntPtr zi_ptr, double b, int m, double eps, IntPtr f_x_ptr)
    {
        double* zr = (double*)zr_ptr.ToPointer();
        double* zi = (double*)zi_ptr.ToPointer();

        f_xa = Marshal.GetDelegateForFunctionPointer<delegatefunc_xa>(f_x_ptr);

        complex z = new complex(*zr, *zi);
        complex r = mtcl(z, b, m, eps);
        *zr = r.R;
        *zi = r.I;
        return 1;
    }

    /// <summary>
    /// 蒙特卡罗法求复根
    /// f计算方程左端函数值的函数名。
    /// </summary>
    /// <param name="z">z根的初值</param>
    /// <param name="b">b均匀分布随机数的端点初值。</param>
    /// <param name="m">控制调节b的参数。</param>
    /// <param name="eps">控制精度要求。</param>
    /// <returns>函数返回根的终值。若程序显示“b调整了100次！迭代不收敛！迭代不收敛！”，则需调整b和m的值再试。</returns>
    public static unsafe complex mtcl(complex z, double b, int m, double eps)
    {
        int flag, k;
        double x1, y1;
        complex z1, zz, zz1;

        k = 0;
        flag = 0;

        //zz = (*f)(z);
        double* zv = stackalloc double[2] { z.R, z.I };
        f_xa(zv);
        zz = new complex(zv[0], zv[1]);

        double* z1v = stackalloc double[2] { 0.0, 0.0 };

        while (flag <= 100)
        {
            k = k + 1;
            // 产生随机数对
            x1 = -b + 2.0 * b * (rnd.NextDouble());
            y1 = -b + 2.0 * b * (rnd.NextDouble());
            z1 = z + new complex(x1, y1);

            //zz1 = (*f)(z1);
            *z1v = z1.R;
            *(z1v + 1) = z1.I;
            f_xa(z1v);
            zz1 = new complex(z1v[0], z1v[1]);

            if (zz1.cfabs() >= zz.cfabs())
            {
                if (k == m)
                {
                    k = 0;
                    flag = flag + 1;
                    b = b / 2.0;
                }
            }
            else
            {
                k = 0;
                z = z1;
                zz = zz1;
                if (zz.cfabs() < eps)
                {
                    return (z);
                }
            }
        }
        return (z);
    }

    /*
    // 实函数方程求复根例
      int main()
      { 
          int m;
          double b,eps;
          complex z, x, mtclf(complex);
          x=complex(0.5,0.5); 
          b=1.0; m=10; eps=0.00001;
          z = mtcl(x,b,m,eps,mtclf);
          cout <<"z = "; z.prt(); cout <<endl;
          cout <<"检验: f(z) = "; mtclf(z).prt(); cout <<endl;
          return 0;
      }
    // 实函数方程
             //要将系数化为复数
    complex mtclf(complex x)
      { 
          complex z;
          z = x*x - complex(6.0,0.0)*x + complex(13.0,0.0);
          return(z);
      }
    */
    /*
    // 复函数方程求复根例
      int main()
      { 
          int m;
          double b,eps;
          complex z, x, mtclf(complex);
          x=complex(0.5,0.5); 
          b=1.0; m=10; eps=0.00001;
          z = mtcl(x,b,m,eps,mtclf);
          cout <<"z = "; z.prt(); cout <<endl;
          cout <<"检验: f(z) = "; mtclf(z).prt(); cout <<endl;
          return 0;
      }
    // 复函数方程
      complex mtclf(complex x)
      { 
          complex z;
          z = x*x + x*complex(1.0,1.0) + complex(-2.0,2.0);
          return(z);
      }
    */
}

